home *** CD-ROM | disk | FTP | other *** search
/ Best of www.BestZips.com (Collector's Edition) / Best of WWW.BESTZIPS.COM Collector's Edition (JCSM Shareware) (JCS Marketing).ISO / prgtools / euphor14.zip / MSET.EX < prev    next >
Text File  |  1996-06-29  |  10KB  |  379 lines

  1.         ---------------------------------
  2.         -- Plotting the Mandelbrot Set --
  3.         ---------------------------------
  4. -- Usage: ex mset [filename.dat]
  5. -- e.g. ex mset msetb
  6.  
  7. -- Either generate the initial picture, or redisplay an old one.
  8. -- Hit Enter at any time to stop the display and then hit Enter again 
  9. -- to display a grid. Use the arrow keys to select the most interesting 
  10. -- box in the grid. Hit Enter to enlarge this box to the full size of 
  11. -- the screen, or hit q to quit. The pictures that you display are saved in
  12. -- mseta.dat, msetb.dat, ... You can redisplay them and then put up the 
  13. -- grid and continue zooming in. As you zoom in, black areas are eroded 
  14. -- around the edges, as the iteration count is increased, and we find that
  15. -- points originally believed to be members of the Mset are not members 
  16. -- after all.
  17.  
  18. without type_check
  19.  
  20. include graphics.e
  21. include select.e
  22. include get.e
  23.  
  24. constant GRAPHICS_MODE = 18  -- Try to increase the graphics mode. 
  25.                  -- See euphoria\include\graphics.e
  26.  
  27. constant ZOOM_FACTOR = 20    -- grid size for zooming in
  28.  
  29. constant FALSE = 0, TRUE = 1
  30. constant REAL = 1, IMAG = 2
  31.  
  32. constant ARROW_LEFT  = 331,
  33.      ARROW_RIGHT = 333,
  34.      ARROW_UP    = 328,
  35.      ARROW_DOWN  = 336
  36.  
  37.     -- types --
  38.  
  39. type natural(integer x)
  40.     return x >= 0
  41. end type
  42.  
  43. type complex(sequence x)
  44.     return length(x) = 2 and atom(x[1]) and atom(x[2])
  45. end type
  46.  
  47. procedure beep()
  48. -- make a beep sound
  49.     atom t
  50.  
  51.     t = time()
  52.     sound(500)
  53.     while time() < t + .2 do
  54.     end while
  55.     sound(0)
  56. end procedure
  57.  
  58. natural ncolors
  59. integer max_color, min_color
  60.  
  61. object prev_mixture 
  62. prev_mixture = {0, 0, 0} -- color 0: black
  63.  
  64. procedure randomize_palette()
  65. -- choose random color mixtures    
  66.     for i = max_color to min_color by -1 do
  67.     prev_mixture = palette(i, rand(repeat(63, 3)))
  68.     end for
  69.     prev_mixture = palette(0, rand(repeat(63, 3)))
  70. end procedure
  71.  
  72. procedure rotate_palette()
  73. -- rotate the color mixtures of all the colors
  74.     
  75.     for i = max_color to min_color by -1 do
  76.     prev_mixture = palette(i, prev_mixture) 
  77.     if atom(prev_mixture) then
  78.         return -- didn't work
  79.     end if
  80.     end for
  81.     prev_mixture = palette(0, prev_mixture)
  82. end procedure
  83.  
  84. natural max_iter
  85.  
  86. sequence vc -- current video configuration
  87.  
  88. procedure grid(sequence x, sequence y, natural color)
  89. -- draw the grid
  90.     atom dx, dy
  91.  
  92.     dx = vc[VC_XPIXELS]/ZOOM_FACTOR
  93.     dy = vc[VC_YPIXELS]/ZOOM_FACTOR
  94.  
  95.     for i = x[1] to x[2] do
  96.     draw_line(color, {{i*dx, y[1]*dy}, {i*dx, y[2]*dy}})
  97.     end for
  98.     for i = y[1] to y[2] do
  99.     draw_line(color, {{x[1]*dx, i*dy}, {x[2]*dx, i*dy}})
  100.     end for
  101. end procedure
  102.  
  103. function zoom()
  104. -- select place to zoom in on next time
  105.     integer key
  106.     sequence box
  107.  
  108.     grid({0, ZOOM_FACTOR}, {0, ZOOM_FACTOR}, 7)
  109.     box = {0, ZOOM_FACTOR-1}
  110.     while TRUE do
  111.     grid({box[1], box[1]+1}, {box[2], box[2]+1}, 14)
  112.     key = get_key()
  113.     if key != -1 then
  114.         grid({box[1], box[1]+1}, {box[2], box[2]+1}, 7)
  115.         if key = ARROW_UP then
  116.         if box[2] > 0 then
  117.             box[2] = box[2]-1
  118.         end if
  119.         elsif key = ARROW_DOWN then
  120.         if box[2] < ZOOM_FACTOR-1 then
  121.             box[2] = box[2]+1
  122.         end if
  123.         elsif key = ARROW_RIGHT then
  124.         if box[1] < ZOOM_FACTOR-1 then
  125.             box[1] = box[1]+1
  126.         end if
  127.         elsif key = ARROW_LEFT then
  128.         if box[1] > 0 then
  129.             box[1] = box[1]-1
  130.         end if
  131.         elsif key >= 27  then
  132.         return {}  -- quit
  133.         else
  134.         return {box[1], ZOOM_FACTOR - 1  - box[2]}
  135.         end if
  136.     end if
  137.     end while
  138. end function
  139.  
  140. procedure save_points(integer file, integer rep_count, integer color)
  141. -- We do a bit of image compression here by recording the number of
  142. -- consecutive points of a certain color. Can have up to 256 colors.
  143.     while rep_count > 255 do
  144.     puts(file, 255)
  145.     puts(file, color)
  146.     rep_count = rep_count - 255
  147.     end while
  148.     if rep_count > 0 then
  149.     puts(file, rep_count)
  150.     puts(file, color)
  151.     end if
  152. end procedure
  153.  
  154. function mset(complex lower_left,  -- lower left corner
  155.           complex upper_right) -- upper right corner
  156. -- Plot the Mandelbrot set over some region.
  157. -- The Mandelbrot set is defined by the equation: z = z * z + C
  158. -- where z and C are complex numbers. The starting point for z is 0.
  159. -- If, for a given value of C, z approaches infinity, C is considered to
  160. -- *not* be a member of the set. It can be shown that if the absolute value
  161. -- of z ever becomes greater than 2, then the value of z will increase
  162. -- towards infinity from then on. After a large number of iterations, if
  163. -- the absolute value of z is still less than 2 then we assume with high
  164. -- probability that C is a member of the Mset and this program will show
  165. -- that point in black.
  166.     complex c
  167.     atom zr, zi, zr2, zi2, cr, ci, xsize, ysize
  168.     natural member, stop, color, rep_count, width, height
  169.     natural file_no
  170.     integer pic_save, prev_color
  171.  
  172.     clear_screen()
  173.     height = vc[VC_YPIXELS]
  174.     width = vc[VC_XPIXELS]
  175.     ncolors = vc[VC_NCOLORS]
  176.     xsize = (upper_right[REAL] - lower_left[REAL])/(width - 1)
  177.     ysize = (upper_right[IMAG] - lower_left[IMAG])/(height - 1)
  178.     c = {0, 0}
  179.  
  180.     -- choose a new file to save the picture into
  181.     file_no = 0
  182.     for i = 'a' to 'z' do
  183.     pic_save = open("mset" & i & ".dat", "rb")
  184.     if pic_save = -1 then
  185.         file_no = i
  186.         exit
  187.     else
  188.         -- file exists
  189.         close(pic_save)
  190.     end if
  191.     end for
  192.     if file_no then
  193.     pic_save = open("mset" & file_no & ".dat", "wb")
  194.     else
  195.     puts(1, "Couldn't find a new file name to use\n")
  196.     return 1
  197.     end if
  198.  
  199.     -- save graphics mode and max_iter
  200.     printf(pic_save, "%d ", vc[VC_MODE])
  201.     printf(pic_save, "%d ", max_iter)
  202.     -- save lower_left & upper_right as floating-point sequences
  203.     printf(pic_save, "{%20.15f,%20.15f}", lower_left)
  204.     printf(pic_save, "{%20.15f,%20.15f}", upper_right)
  205.     max_color = -1
  206.     min_color = 99999
  207.     for y = 0 to height - 1 do
  208.     if get_key() != -1 then
  209.         close(pic_save)
  210.         return 0
  211.     end if
  212.     c[IMAG] = upper_right[IMAG] - y * ysize
  213.     prev_color = -1 -- start fresh for each line
  214.     rep_count = 0
  215.     for x = 0 to width - 1 do
  216.         c[REAL] = lower_left[REAL] + x * xsize
  217.         member = TRUE
  218.         zr = 0
  219.         zi = 0
  220.         zr2 = 0
  221.         zi2 = 0
  222.         cr = c[REAL]
  223.         ci = c[IMAG]
  224.         for i = 1 to max_iter do
  225.         zi = 2 * zr * zi + ci
  226.         zr = zr2 - zi2 + cr
  227.         zr2 = zr * zr
  228.         zi2 = zi * zi
  229.         if zr2 + zi2 > 4 then
  230.             member = FALSE
  231.             stop = i
  232.             exit
  233.         end if
  234.         end for
  235.         if member = TRUE then
  236.         color = 0
  237.         else
  238.         color = stop + 51 -- gives nice sequence of colors
  239.         while color >= ncolors do
  240.             color = color - ncolors
  241.         end while
  242.         if color > max_color then
  243.             max_color = color
  244.         end if
  245.         if color < min_color then
  246.             min_color = color
  247.         end if
  248.         end if
  249.         pixel(color, {x, y}) -- also show non-member "bands"
  250.         if color = prev_color then
  251.         rep_count = rep_count + 1
  252.         else
  253.         save_points(pic_save, rep_count, prev_color)
  254.         rep_count = 1
  255.         prev_color = color
  256.         end if
  257.     end for
  258.     -- close off count at end of each line
  259.     save_points(pic_save, rep_count, color)
  260.     end for
  261.     beep()
  262.     close(pic_save)
  263.     return 0
  264. end function
  265.  
  266. procedure view(integer pic_save)
  267. -- redisplay a saved picture file
  268.     integer x, color, rep_count
  269.  
  270.     max_color = -1
  271.     min_color = 99999
  272.     for y = 0 to vc[VC_YPIXELS] - 1 do
  273.     x = 0
  274.     while x < vc[VC_XPIXELS] do
  275.         rep_count = getc(pic_save)
  276.         color = getc(pic_save)
  277.         if rep_count <= 0 then
  278.         return
  279.         end if
  280.         if rep_count = 1 then
  281.         pixel(color, {x, y})  -- faster
  282.         x = x + 1
  283.         else
  284.         draw_line(color, {{x, y}, {x+rep_count-1, y}})
  285.         x = x + rep_count
  286.         end if
  287.         if color != 0 then
  288.         if color > max_color then
  289.             max_color = color
  290.         end if
  291.         if color < min_color then
  292.             min_color = color
  293.         end if
  294.         end if
  295.     end while
  296.     end for
  297. end procedure
  298.  
  299. procedure Mandelbrot()
  300. -- main procedure
  301.     sequence delta, new_box
  302.     complex lower_left, upper_right
  303.     sequence cl, dataname, g
  304.     integer pic_save, mode
  305.  
  306.     cl = command_line()
  307.     if length(cl) >= 3 then
  308.     -- redisplay a saved picture
  309.     dataname = cl[3]
  310.     pic_save = open(dataname, "rb")
  311.     if pic_save = -1 then
  312.         if not find('.', dataname) then
  313.         dataname = dataname & ".dat"
  314.         pic_save = open(dataname, "rb")
  315.         end if
  316.         if pic_save = -1 then
  317.         puts(1, "Couldn't open " & dataname & '\n')
  318.         return
  319.         end if
  320.     end if
  321.     g = {}
  322.     for i = 1 to 4 do
  323.         g = g & get(pic_save)
  324.     end for
  325.     if g[1] != GET_SUCCESS or
  326.        g[3] != GET_SUCCESS or
  327.        g[5] != GET_SUCCESS or
  328.        g[7] != GET_SUCCESS then
  329.         puts(1, "Couldn't read " & dataname & '\n')
  330.         return
  331.     end if
  332.     mode = g[2]
  333.     max_iter = g[4]
  334.     lower_left = g[6]
  335.     upper_right = g[8]
  336.     if graphics_mode(mode) then
  337.     end if
  338.     vc = video_config()
  339.     view(pic_save)
  340.     else
  341.     -- initially show the upper half:
  342.     max_iter = 30 -- increases as we zoom in
  343.     lower_left = {-1, 0}
  344.     upper_right = {1, 1}
  345.     -- set up for desired graphics mode
  346.     if not select_mode(GRAPHICS_MODE) then
  347.         puts(2, "couldn't find a good graphics mode\n")
  348.         return
  349.     end if
  350.     vc = video_config()
  351.     if mset(lower_left, upper_right) then
  352.         return
  353.     end if
  354.     end if
  355.  
  356.     while TRUE do
  357.     while get_key() = -1 do
  358.         --rotate_palette()
  359.         randomize_palette()
  360.     end while
  361.     new_box = zoom()
  362.     if length(new_box) = 0 then
  363.         exit
  364.     end if
  365.     delta = (upper_right - lower_left)
  366.     lower_left = lower_left + new_box / ZOOM_FACTOR * delta
  367.     upper_right = lower_left + delta / ZOOM_FACTOR
  368.     max_iter = max_iter * 2  -- need more iterations as we zoom in
  369.     if mset(lower_left,  upper_right) then
  370.         exit
  371.     end if
  372.     end while
  373. end procedure
  374.  
  375. Mandelbrot()
  376.  
  377. if graphics_mode(-1) then
  378. end if
  379.